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ABSTRACT 

We study linear and nonlinear development of relativistic and ultrarelativistic 
current sheets of pair (e^) plasmas with antiparallel magnetic fields. Two types of 
two-dimensional problems are investigated by particle-in-cell simulations. First, 
we present the development of relativistic magnetic reconnection, whose outflow 
speed is on the of the light speed c. It is demonstrated that particles are strongly 
accelerated in and around the reconnection region and that most of the magnetic 
energy is converted into a "nonthermal" part of plasma kinetic energy. Second, 
we present another two-dimensional problem of a current sheet in a cross-field 
plane. In this case, the relativistic drift kink instability (RDKI) occurs. Parti- 
cle acceleration also takes place, but the RDKI quickly dissipates the magnetic 
energy into plasma heat. We discuss the mechanism of particle acceleration and 
the theory of the RDKI in detail. It is important that properties of these two 
processes are similar in the relativistic regime of T > mc 2 , as long as we consider 
the kinetics. Comparison of the two processes indicates that magnetic dissipa- 
tion by the RDKI is a more favorable process in the relativistic current sheet. 
Therefore, the striped pulsar wind scenario should be reconsidered by the RDKI. 

Subject headings: acceleration of particles — magnetic fields — plasmas — rela- 
tivity — stars: winds, outflows — pulsars: individual (Crab Pulsar) 

1. Introduction 

Magnetic reconnection is an important large-scale processes in space plasmas. By rear- 
ranging the magnetic field topology, it causes heating and particle acceleration of plasmas as 
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well as dissipation of the magnetic fields. Although it has been extensively studied in wide 
range of solar terrestrial sites — stellar and solar flares (Parker 1957), the planetary mag- 
netospheres (Dungey 1961 Russell et al. 1998) and solar wind (Phan et al. 2006) — it has 



also been discussed in high-energy astrophysical contexts such as the magnetized loop of the 



Galactic center 


(Heyvaerts et al. 1988), the jets from active galactic nuclei (di Matteo 


1998 


Schopper et al. 


1998 Lesch & Birk 


1998 Larrabee et al. 2003 


, quite probably gamma-ray 


bursts (Drenkhahn & Spruit 2002) 


, and pulsar winds (Michel 


1982 Coroniti 1990; 1 


vlichel 



1994). In these situations, reconnection often takes place in relativistic electron-positron 
plasmas. Especially in the Crab pulsar wind, in the relativistic outflow of pair plasmas from 
the central neutron star, relativistic reconnection or relevant processes have been considered 
as possible processes to explain the long-standing "er-problem" (<r is the ratio of the Poynting 
flux energy to the particle kinetic flux); the plasma outflow is originally Poynting-dominated 
(a ~ 10 4 ) close to the neutron star (Arons 1979), but it should be kinetic-dominated near 
the downstream termination shock 



[a 
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Kennel & Coroniti (1984)). The possible dis 



sipation mechanisms are instabilities in the striped current sheets; due to the fast rotation 
of the central neutron star, whose magnetic dipole is in the oblique direction, alternating 
magnetic fields are strongly "striped" near the equatorial plane. We assume that there are 
current sheets between such opposite magnetic field lines, and magnetic reconnection occurs 
and dissipates the magnetic energy there (Coroniti 1990). Based on this idea, magnetic 



dissipation in the pulsar wind has been estimated by using a one-dimensional MHD model 
(Lyubarsky & Kirk 2001 Kirk & Skjaeraasen 2003), but they were not so successful. The 



main problem is that basic properties of relativistic reconnection or current sheet processes 
are still unclear. 

Until recent years, there have been few theoretical studies on relativistic reconnection. 



Zelenyi & Krasnoselskikh (1979) studied the relativistic kinetic description of the tearing 



instability, which is the most important instability in the reconnection context. Blackman Sz 



Field (1994) studied the steady state reconnection models in relativistic pair plasmas. They 
argued that the reconnection rate becomes closer to the unity due to Lorentz effects, and 
therefore, faster energy conversion is possible. These steady models are further investigated 



by several authors (Lyutikov & Uzdensky 2003 Lyubarsky 2005), although there is contro- 



versy about whether inflow speed can be relativistic or not. On the viewpoint of particle 



acceleration, Romanova & Lovelace (1992) studied the particle orbits in the reconnection 
field model, and they obtained an energetic energy spectrum of pair plasma, which is ap- 
proximated by the power-law distribution with an index of —1.5. Larrabee et al. (2003) 
obtained the power-law index of —1, assuming that particle motion is restricted in the neu- 
tral plane. In a fully self-consistent way, relativistic magnetic reconnection has recently been 
explored by particle-in-cell (PIC) simulations ( Zenitani fc Hoshino||2001 2005b Jaroschek et 



al. 2004). Zenitani & Hoshino (2001) demonstrated that an enhanced dc acceleration takes 
place in and around the X-type region, due to the reconnection electric field E y . This accel- 
eration is so strong that a power law index of the energy spectrum around the acceleration 
site is on the order of —1 (Zenitani & Hoshino 2001), and the energy spectrum over the 



whole simulation domain is approximated by the power-law with an index of —3 (Jaroschek 



et al.|200~4~ ). Using relativistic resistive MHD code, Watanabe & Yokoyama (2006) presented 
Petschek reconnection in the mildly relativistic case. 

The current sheet configuration with antiparallel magnetic field lines is also influenced 
by cross field instabilities, whose wavevectors are perpendicular to the reconnection plane. 
In geophysical contexts, the lower hybrid drift instability (LHDI) ( |Krall fc Liewer 1971 



Davidson et al. 1977) may be of importance, and it leads to transport of plasma heat (Huba 



& Drake 1981) and the fast triggering of magnetic reconnection (Shinohara & Fujimoto 



2005 Tanaka et al. 2006). The Kelvin-Helmholtz instability (KHI) (Yoon fc Drake 1996 



Shinohara et al.|2001 ), which arises from the velocity shear between the fast-drifting plasmas 
in the current sheet and the background plasmas, is of importance because of its significant 
modulation, and it enhances the magnetic diffusion rate by the LHDI ( Shinohara et al.|2001 ). 
Recently the drift kink instability (DKI) ( |Zhu & Winglee||1996| |Pritchett et"aT|l996p and its 
cousin mode of the drift sausage instability (DSI) dBiichner & Kuika1[l999| |Yoon & LuipOOl 



Silin &: Buchner||2003 ) have been introduced. The DKI is a long- wavelength electromagnetic 
mode, which is driven by the fast-drifting plasmas in a thin current sheet and it quickly grows 



when the mass ratio of the positively charged particles to the electrons is small (Daughton 



1998, 1999). Therefore the DKI smears out in geophysical ion-electron situations but the 



nonlinear development of the LHDI may lead to the current sheet modulation by the DKI 



(Horiuchi & Sato 1999). On the contrary, the DKI quickly grows in relativistic pair plasmas 



in which the mass ratio is the unity, while the LHDI smears out. Zenitani & Hoshino (2003 



2005a) investigated the DKI in relativistic current sheets of pair plasmas, and they found a 



dc particle acceleration along the neutral sheet in the nonlinear stage of the relativistic drift 
kink instability (RDKI). Because of its fast growth rate, the RDKI is one of the most likely 
processes in relativistic current sheets. 

The purpose of this paper is to investigate basic physical properties of reconnection or 
reconnection-related processes in relativistic pair plasmas. Two types of two-dimensional 
problems, relativistic magnetic reconnection and the RDKI, are demonstrated by using fully 
electromagnetic PIC simulations. In addition, based on the obtained theories we discuss 
possible application to the Crab pulsar wind problem. 

This paper is organized as follows. In §2 we describe our simulation setup. In §3 a 
two-dimensional simulation study of relativistic magnetic reconnection is presented. In §4 
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another two-dimensional study of the RDKI is presented. In §5 we compare the two results 
and discuss the possible application to the pulsar winds. 



2. Simulation 

2.1. Simulation method 

The simulation is carried out by a three-dimentional electromagnetic PIC code. In this 
code we solve the following basic equations: the relativistic equation of motion 

= 9j (e + ^xb) (1) 

d . . 

- Xj = Vj (2) 



dt 3 



and the Maxwell equations 



VXB - £ J+ i». ,3) 

1 9 

VxE = — -B (4) 

C OT 

V-B = (5) 
V-E = 4vrp (6) 

where Xj, v,, 7^, m^, and ^ are the position, the velocity, the Lorentz factor, the rest mass, 
and the charge for the j-th super-particle in the simulation system. The charge density p 
and the current density j are obtained by 

3 
3 

where denotes integration in the cell and S(x) is the shape function of the super-particles. 
For simplicity, we do not consider any collisions, radiation, pair production, or pair annihi- 
lation of electron-positron plasmas. 



2.2. Initial configuration 



As an initial configuration, we use a relativistic Harris configuration (Kirk & Skjaeraasen 
2003 Harris 1962). In the Cartesian coordinate system, magnetic field and plasma distribu- 
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tion functions are 



B = £ tanh(z/A)x, (9) 

n cosh" 2 (z/ A) -7 s (e - f3 s mcu y ) , 

h ~ Airm 2 cTK 2 (mc 2 /T) eXlp[ T 1 

exp[-^], (10) 



ATcm 2 cT bg K 2 (mc 2 /T b g) T ht 

where A is the thickness of the current sheet, the subscript s denotes the species (plus signs 
for positrons and minus signs for electrons) , no is the plasma number density of the current 
sheet in the proper frame, T is the plasma temperature including the Boltzmann constant, 
K 2 (x) is the modified Bessel function of the second kind, cj3 s is the drift speed of the species, 
and u is the four- velocity. Throughout this paper we set (3 + = 0.3 and /?_ = —0.3 and 
therefore 7/3 = 1.048. The is the number density of background plasmas, and Tb g is 
its temperature including the Boltzmann constant. The pressure balance condition and the 
current condition are satisfied in the equilibrium state; 

B 2 /8n = 2d T/^ = 2n T (11) 

cB / (4ttA) = ^ lsq s n s v s = 2en^ p cf3 = 2ed cp, (12) 

s 

where do = 7,3/10 is the plasma density in the laboratory frame. 

In the Harris system, we can take two free parameters. We employ the temperature 
T and the drift speed parameter f3. The temperature T is a measure of the relativity in 
this study. In the nonrelativistic regime, one can obtain the typical Alfven speed Va in the 



system from equation (11). 



y / 47rm(2n ) " " Vmc 2 . 

Therefore, the typical Alfven speed becomes on the order of c and then several relativistic 
effects appear, when T > mc 2 . The drift speed parameter (3 also stands for the current sheet 



thickness. Combining equations (11) and (12), one can obtain the Debye length Ad 



^'i^r^- < 14 > 

The typical gyroradius tl in the T 3> mc 2 limit can be approximated as 

feBoy 1 (T 7lf3)mc 2 

r L ~ c = — = 03/27/3 • A. 15 

K^mcJ mc 1 eB 
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In the nonrelativistic limit of T mc 2 , 

f eB yi mc /~2T~ (3 / T y 1 ' 2 x ,, R s 
V mc / ei3 V m 7/3 J2 r )^\mc A / 

where is the thermal velocity of plasmas. 

We use the three-dimensional PIC code. The system size is set to 1024 (x) x 1 (y) x 512 
(2) in reconnection studies and 1 (x) x 256 (y) x 512 (z) in drift kink studies. We consider 
periodic boundaries in the x-, y- and z-directions. Since the magnetic field lines change their 
directions in a current sheet, we set two simulation domains in the z direction: one domain 
for the first Harris sheet and the other for the second Harris sheet, which has the same 
physical properties in the opposite directions. Usually, physical processes are investigated in 
the first half (bottom half) of the whole simulation box, and so the effective size in z is 256 
grids. The typical scale of the current sheet A is set to 10 grids. Therefore, the boundaries 
are located at x — ±51. 2A (reconnection) or at y — ±12. 8A (drift kink). The ^-boundaries 
of the main simulation domain are located at z = ±12. 8A. Simulation time is normalized by 
the light transit time r c = A/c. 

A list of simulation runs is presented in Table [3] Note that tu c = cBq / mc means "unit 
gyroradius". In the relativistic limit, the typical gyroradius becomes larger by a factor of T 
(eq. (15)). 

3. Relativistic magnetic reconnection 

3.1. Simulation result 

First, we present our simulation study of relativistic magnetic reconnection. The system 
size is 1024 (x) X 1 (y) X 512 (z), and the physical size of the main simulation box is —51.2 < 
x/X < 51.2, —12.8 < z/X < 12.8. The plasma temperature is set to T/mc 2 = 1.0, so that 
the system's typical Alfven speed is Va ~ 0.53c. The background plasma temperature T& g is 
set to O.lmc 2 . We call this run run R3 in Table [3] In order to trigger an X-type neutral line 
around the center of the simulation box, we add artificial electric fields in the very early stage 
of the simulation. We assume the triggering electric field E tr i g and its maximum amplitude, 
typical location, and duration time are 0.3B q (Va/c), (x ± Ax, z ± Az) = (0 ± 2, 3 ± 1), and 
10A/Va, respectively. During < t < T trig) we force E y > E tr i g so that plasmas enter the X 
point. Its duration 10A/Va ~ 10 — 15r c is short enough, compared with the total timescale 
of the simulation. Because of this modification the system slightly gains energy by 0.2% in 
the very early stage, and the total energy is conserved within 0.4% error until the end of the 
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simulation at t/r c = 300.0. 

The three panels in Figure [l] present the snapshots at t/r c = 60.0,80.0, and 100.0. 
The color contour shows the density of plasmas, which is normalized by the original plasma 
sheet density po — 27/^0. The solid lines represent magnetic field lines. The field lines are 
originally in the +x direction on the upside of the current sheet and the — x direction on the 
downside of the current sheet. Because of the triggering field, the field lines start to reconnect 
near the center of the simulation box. Along with the magnetic field lines, plasmas from the 
background region come into the X-type region and then they stream out as reconnection jets 
toward the ±x directions. The evolution looks like Sweet-Parker type reconnection with a flat 
current sheet structure. The average velocity of plasmas (v) = J /vrfu/ J fdu is up to (0.8- 
0.9)c in the outflow regions and (0.3-0.4)c in the inflow regions. The outflow speed exceeds 
the typical Alfven speed in the system. At t/r c = 80.0 and 100.0, reconnected magnetic 
field lines {B z ) are swept away from the A-type region. As a result, the field lines are piled 
up around x/X = ±20 or ±30 in front of the dense plasma regions. After t/r c = 100.0, 
the system is influenced by the periodic boundary conditions. Reconnection outflow jets 
come closer to the x-periodic boundaries at x/X = ±51.2, where the plasma density starts to 
increase. In addition, the A-point starts to sweep magnetic flux from the other simulation 
box. At t/r c = 100.0, the magnetic field lines at the A-point come from z/X ~ ±12-13 
around x/X = ±51.2. Note that the z-boundaries are located at z/X — ±12.8. Dense 
plasma regions are blown away from the center, they meet each other near the x-periodic 
boundaries. Then they merge into a single O-point, and the O-point starts to evolve in 
the vertical (z) direction. After that, the system evolves gradually to another relaxed state: 
vertical current sheets. Figure [2] is a snapshot in the late-time stage at t/r c = 300.0. Both 
main (—12.8 < z/X < 12.8) and sub (12.8 < z/X < 38.4) simulation boxes are presented. 
The white rectangle indicates the region in the snapshots from Figure [TJ The magnetic fields 
are antiparallel; B z < in the left half (z/X < 0) and B z > in the right half (0 < z/X), but 
at present we do not observe the secondary magnetic reconnection in these vertical current 
sheets. 

The panels in Figure [3] are energy spectra in the main simulation box. The vertical 
axis shows the count number of super particles on a log scale, and the horizontal axis shows 
particle's energy which is normalized by the rest-mass energy mc 2 . The bottom panel of 
Figure [3] presents the energy spectra in double-log format. The spectrum labeled by t/r c = 
100.0' is the partial energy spectrum at t/r c = 100.0, which is integrated near the central 
region of —12.8 < x/X < 12.8. After the reconnection breaks up, the high-energy tail rapidly 



grows in time. This is due to the dc particle acceleration around the A-point (Zenitani & 



Hoshino 2001) and we discuss the acceleration mechanism in more detail in §3.2. The slope 



of the spectrum at t/r c = 80.0 is well approximated by the power-law distribution with an 



-8- 



index of —3.1 to —3.2 in the range of 10 < e/mc 2 < 50. This is consistent with the other 
study; Jaroschek et al. (2004) obtained a spectral index of —3 in their PIC simulations with 
lower temperature T/mc 2 = 0.2. At t/r c = 100.0, the spectrum seems to have a double 
power-law shape, whose hardest slope is roughly approximated by the index of —2.4. The 
partial energy spectrum at t/r c = 100.0 looks harder. Finally, the maximum energy is up to 
150mc 2 . The late-time state at t/r c = 300.0 also has its spectral index of —2.4. 



The panels in Figure [4] show the electric field structure at t/r c = 80.0 and 100.0. Because 
of the symmetry, we present two properties, (E 2 — B 2 ) in the left halves (x/X < 0) and the 
reconnection electric field E y in the right halves (0 < x/X). The white contour lines are 
drawn by 0.25. Since the reconnection magnetic fields (B x or B z ) are small near the X-type 
region and since the reconnection electric field E y is finite, (E 2 — B 2 ) becomes positive inside 
the narrow sheet region around the X-type region; —10 < x/X < 10 at t/r c = 80.0 and 
—15 < x/X < 15 at t/r c = 100.0. The reconnection electric field E y is almost constant 
around the X-type region: E y /B ~ 0.2 at t/r c = 80.0 and E y /B ~ 0.15 at t/r c = 100.0. 
There are strong peaks around x/X = ±20 at t/r c = 80.0, x/X = ±30 at t/r c = 100.0. These 
regions are identical to the magnetic pile-up regions, where the reconnected magnetic field 
lines are compressed in front of the dense current sheet regions. The electric fields E y are 
enhanced due to the motional electric field of the magnetic piled-up flux. The other two 
components (E x and E z ) are negligible. 



3.2. Particle acceleration 

Next, we look at the particle acceleration in relativistic reconnection. In order to analyze 
how and where the nonthermal particles are accelerated, we pick up high-energy particles 
whose energy exceeds 50mc 2 at t/r c = 100.0. The panels in Figure [H] show their spatial 
distribution at three different stages. At t/r c = 40.0, some particles are near the X-type 
region and other particles are widely scattered in the background region. Importantly, they 
stay around or move into the X-type region at t/r c = 60.0, when the reconnection breaks up. 
The nonthermal tail in the energy spectra starts to grow in this stage. As the reconnection 
evolves, they are spread over the thin current sheet, and then some of them are ejected from 
the reconnection region. The acceleration site is obviously the central region in and around 
the X-type region. Looking at the electromagnetic fields near the acceleration site in panels 
in Figure [4], we note that the electric fields are relatively stronger than the magnetic fields. 

We call the central region that satisfies (E 2 — B 2 ) > the "acceleration region" (AR), 
where the electric field is dominant. In two-dimensional reconnection, (E 2 — B 2 ) > is 



equivalent to our old definition of the AR: \E y \ > \B 2 ± B 2 \ in Zenitani & Hoshino (2001) 
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(hereafter ZH01). Such an AR does not have a de Hoffman- Teller frame. For example, let us 
consider an electromagnetic field B(0, 0, B x ), E(0, E y , 0) near the X-points along the neutral 
sheet (z = 0). The fields can be transformed into 

R(nnB , pm pnW | B'(0,0,B,/ 7w ),E'(0,0,0) (\E y \ < \B,\) 
B(0,0,BJ,E(0,B s ,0) | B , (0i0i0)iE , (0 _ v ^ |0) (|Bs| > m (17) 

where the prime is for the transformed properties, v\/c = E y /B z and vi/c = B z /E y are 
the transformation speeds, and j Vl and'y V2 are the relevant Lorentz factors. Generally, when 
the Lorentz invariant (E 2 — B 2 ) is positive, the electric field cannot be removed and the 
particle motion in the AR is controlled by the electric field in the transformed frame. In 
the observer's frame, particles are driven into some ^/-aligned direction, gaining their energy 
through the reconnection electric field E y . So, once particles enter the acceleration site, 
they are driven by the electric field E' y inside the AR and they travel through (relativistic) 
Speiser orbits (Speiser 1965) outside the AR. The AR looks thinner in z compared with 
ZH01, but the thickness of the AR does not provide significant change because the magnetic 
field B x confines particles near the neutral sheet in a reconnection configuration. In ZH01, 
we estimate that the energy distribution in the AR can be approximated by the power law 
distribution, whose index is on the order of —1. Although it is difficult to evaluate the partial 
energy spectra in the AR by power law, the index is between —2 and —1 in Figure [3} 

In the bottom panel of Figure [5j the red spots present 33 particles, whose energy exceeds 
90mc 2 at t/r c = 100.0. They are found only near the magnetic pileup regions along the flat 
current sheet layer. According to the partial and the main energy spectra at t/r c = 100.0 
in Figure [3j the intermediate-energy particles (40mc 2 < e < 70mc 2 ) are found both inside 
and outside the selected region and their spectral slopes look similar. The selected region 
has fewer intermediate-energy particles by a factor of 3, because the region is narrower 
than the typical current sheet length. On the contrary, high-energy particles whose energy 
exceeds e > 80mc 2 are only found outside the selected region. Therefore, these highest 
energy particles are accelerated outside the selected region of —12.8 < x/X < 12.8. In fact, 
during t/r c = 60.0 and 80.0, the highest-energy shoulder increases by Ae ~ (25 — 30)mc 2 , 
and some particles apparently gain more energy than the maximum energy gain inside the 
AR, Ae ~ {eEyC ■ 20r c ) ~ 25mc 2 , where we use E y < 0.2Bq as the typical value. Then, 
we traced the most energetic positron, whose energy is 97mc 2 at t/r c = 100.0 and who is 
found at (x/X,z/X) = (31.0, —0.64). Its trajectory and xz-jxy- projections are presented in 
Figure [6] Labels in the xy plane indicate the relevant times. Time histories of the relevant 
properties are presented in Figure [7j too. Originally, the positron's motion is well explained 
by the meandering motion; it bounces around the neutral sheet in the ±z direction. Around 
t/r c = (40 ~ 50), the reconnection breaks up and the positron enters the A-type region 
(x/X ~ 0). Then, it is driven by the reconnection electric field E y and its energy starts 
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to increase. In the gray background regions, the field satisfies the condition of the AR, 
(E 2 — B 2 ) > 0. This particle does not always stay in the AR, because it travels near the 
edge of the AR. After t/r c = 65, the positron departs from the central AR and then it 
starts to travel through the (relativistic) Speiser orbit. The energy already reaches up to 
e = 50mc 2 through particle acceleration in/near the central AR. At this stage, the positron's 
Larmor radius is ( , ymc 2 /eB ) ~ 8A. Because of a relativistic inertia effect, it can travel a 
longer distance into the y-direction through Speiser orbit. If we observe this particle in the 
moving frame of the plasma outflow, its energy stays nearly constant after t/r c > 65. In 
the simulation frame, the particle gains further energy from E y until it is ejected into the 
x-direction. 

Interestingly, in this case the particle gains more energy than the standard Speiser orbit 
case, because the electric field E y increases in time at the particle's position. Figure [8]is an x- 
t diagram of the positron's position and the .E^-profile. The dashed line presents the positron 
trajectory, and the stack plots show the time development of the reconnection electric field 
E y along the neutral sheet [z = 0). Surprisingly, the positron always stays around the peak 
of the electric field E y , which is related to the magnetic pile-up region. Since magnetic field 
lines are continuously transported from the X-type region, the piled-up peak increases in 
its height, and then E y is further enhanced by the relevant motional electric field. In other 
words, the particle speed resonates with the propagating speed of the reconnection pileup 
front ~ 0.6 — 0.7c, which is slightly slower than the reconnection outflow speed. In the 



geophysical ion-electron contexts, Birn et al. (1997) studied ion acceleration and Hoshino 



et al. (2001) studied electron acceleration in the magnetic dipole/pileup regions. In both 
cases, particles gain their energy from E y , through VB drift motion in the y-direction. In 
our case, the particle directly resonates with the fields. In order to stay for a long time in 
the piled-up front, the particle has to be heavy enough or the particle's Larmor radius has 
to be large enough. Therefore, only high-energy particles can be continuously accelerated 
to the highest energy; most of them are preaccelerated around the AR when they enter the 
piled-up regions. Low-energy particles are quickly ejected into the outflow directions. When 
these pileup regions hit the downstream dense regions, it releases a lot of energy into the 
downstream plasma energy. The downstream plasmas increase their energy on the order of 
As ~ 2r' L eE' y ~ 3(7?nc 2 / eBo)eBo ~ 37mc 2 , where the prime denotes the physical properties 
in the piled-up front frame, and this mirror effect explains the mild enhancement of the 
global energy spectra around e ~ (20 — 30)mc 2 . If we compare the piled-up region with the 
AR, the piled-up region mainly contributes to the plasma heating, but it also enhances the 
maximum energy of a small number of highly accelerated particles. A similar mechanism is 



reported in Jaroschek et al. (2004), although their difference between the dc acceleration in 



the centeral AR and this piled-up acceleration was not clear. Their system becomes highly 
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turbulent, because the main reconnection is suppressed by the conductive wall boundary 
effect in the inflow region. On the contrary, in our case, the main reconnection is so powerful 
that it sweeps out everything into the ±x direction, including small tearing islands in its early 
stage. Therefore, the acceleration continues for a longer time, until outflow jets reach the 
periodic boundaries. We do not know how long the single reconnection region can dominate, 
but we remark that the higher maximum energy will be obtained when a single reconnection 
dominates. 

Figure [9] summarizes the acceleration mechanism, presenting the electromagnetic field 
properties at t/r c = 100.0 along the neutral plane z = 0. As stated, the electric field is 
larger than the magnetic field B z near the X-type region x/X < 14, while B z exceeds E y in 
x/X > 14. The dashed line presents the ratio of E y to B z , which agrees with the plasma 
outflow speed around the outflow the region. In this figure, we can classify the acceleration 
site into the following three types: 

Central acceleration region — Where (E 2 — B 2 ) > and dc direct acceleration takes 
place. 

Outside region — Where (E 2 — B 2 ) < and particles are accelerated through the (rela- 
tivistic) Speiser orbit. 

Magnetic pileup region — Where some high-energy particles resonate with the enhanced 
electric field. 

Regarding the boundary effects on the late-time development, the plasma collision across 
the x-boundaries enhances the magnetic pileup flux and the relevant acceleration. However, 
such a collision is quite likely to occur in realistic situations. For example, the coalescence 
of multiple reconnection islands have recently attracted our attention because it generates 



energetic particles (Saito & Sakai 2006). Since our system width (L x = 102. 4A) is rather 
larger than the typical scale of the tearing mode (~ 20A, assuming that the fastest mode 
is k x X ~ 0.3), we believe that the x-boundary effect is plausible. The z-boundary effect 
may be rather artificial; however, it is unlikely that the z-boundary effect accelerates the 
reconnection process in the early stage of t/r c < 100. The inflow magnetic flux continually 
decreases in time (approximately half of the original flux ~ 0.5-Bo at t/r c = 100.0) and then 
it starts to increase in the later stages at t/r c > 170, due to the O-point expansion in the sub 
simulation domain. In the striped pulsar wind condition, it is possible that double/multiple 
current sheets are close to each other. Since we note that the late-time energy spectra become 
even harder (with an index of ~ —2.4), these results imply that harder energy spectra may 
be generated when the pulsar wind is crowded by the multiple current sheets. Of course, it is 
very difficult to discuss the spectral index from the simulation results, because it shows the 
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very early development in a small region compared with actual astrophysical phenomena. 
It just gives us a hint to discuss an astronomical long-time/large-scale evolution and the 
resulting energy spectra. 



3.3. Case studies 



We have further carried out five simulation runs with different temperature parameters; 
T/mc 2 = 16, 4, 1, 1/4 and 1/16. We call them runs R1-R5 and their parameters are presented 
in Table [3] The background plasma temperature T bg is always set to 0.1T and ^3/(7^0), 
the ratio of the number density in the simulation frame, is fixed to 0.05. The amplitude 
and duration of triggering fields are selected case by case, but we usually use smaller values 
than run R3. For comparison, our previous work (ZH01) is presented as run R0, whose 
parameters are equivalent to T/mc 2 ~ 1/4 and (3 ~ 0.27. 



The results of the parameter studies are presented in Figures 10 11, and 12 and Table [TJ 
Figure 10 shows typical snapshots in runs Rl and R5 (T/mc 2 = 16 and 1/16). The energy 
spectra in run Rl are also presented in Figure fTTl The enhancements of plasma kinetic 



energies in runs R1-R5 are presented in Figure 12 s. Their energies are normalized by the 
plasma pressure energy in the original current sheet, which is roughly proportional to the 
temperature T. The time axis is re-arranged so that reconnection breaks up nearly at the 
same time. In addition, in order to evaluate the acceleration amplitude, we calculate the 



"nonthermal ratio" index of K nt h/K in Figure 12 b, where K is the total kinetic energy in 
the main simulation box and K nt h is its nonthermal part. We estimate K nt h by comparing 
the energy spectra and its equivalent "thermal" distribution. Its derivation is described in 



Appendix A. Figure [12] c presents the size of the acceleration region in the main simulation 
box, which satisfies the condition of (E 2 — B 2 ) > 0. Table [l] shows characteristic speeds in 
the system. We pick up the maximum plasma speed along outflow/inflow regions, as long 
as it does not violate the MHD frozen-in condition so much. The typical Alfven velocity 
is calculated from the background magnetic field Bq and plasma parameters in the current 
sheet (no, T). The inflow Alfven velocity is the local Alfven speed in the background region. 

System evolutions are quite similar in relativistic cases of T/mc 2 > 1. In run Rl, the 



spatial structure (Fig. 10 a) and the power-law index of —3.2 (Fig. 11) in the reconnection 
phase are in good agreement with those of run R3. In the very late time stage (at t/r c = 
314.0), the power-law index is approximately —2.4. The evolution timescale of (30— 50)r c and 
the saturation level (Fig. 12 a) and a portion of the nonthermal energy (Fig. [l2~]fr) are almost 



the same, too. The ratio of the inflow speed and the outflow speed usually stay around 0.3- 
0.4, except for run Rl (Table Q. In run Rl, the rate becomes slightly higher due to the faster 



13 



inflow speed. However, it still stays nonrelativistic; the relevant Lorentz factor is only ~ 1.25, 
and the typical inflow speed is slower because we picked up the highest inflow speed along 
the y axis. Considering that the system is magnetically dominated [(Bq/8ti)/(J2 £ ) ~ 60] in 
the inflow region and that the reconnection layer's aspect ratio is L/5 ~ 10, where L and 
5 are the typical width and height of Sweet-Parker reconnection layer, our results do not 



agree with Lyutikov & Uzdensky (2003); they predicted that the inflow velocity becomes 
relativistic and that the outflow velocity is faster than the inflow Alfven speed when the 
inflow region is magnetically dominated. 



In less relativistic cases, the current sheet structure is ambiguous or diffusive (Fig. [I0p). 
This is due to relatively large Larmor radius (eq. Jl6] ; ~ 0.8A in run R5) in the non- 
relativistic limit. In addition, the system evolves slowly because of the slower characteristic 



speed: the nonrelativistic Alfven speed. Contrary to relativistic MHD simulations (Watan- 
abe fc Yokoyama|2006 ) , in which the outflow speed is well-approximated by the inflow Alfven 
speed, the maximum outflow speed stays between the two Alfven speeds. The size of the 
AR decreases and then fewer nonthermal particles are observed. Strictly speaking, this may 
contain ARs in several X-points or magnetic nulls inside the O-points. However, in the 
breakup stage of reconnection (50 < t/r c < 100), they represent the size of the AR in the 
single reconnection region very well. We believe that the size of the AR is controlled by the 
reconnection outflow speed. The magnetic topology of reconnection is rather insensitive to 
parameters and so the relative size of the AR depends on the amplitude of the reconnection 
electric field E y . Then, E y is controlled by the typical outflow speed through the frozen- in 
condition at the outflow region, \v out \ = c\E y \/\B z0 \, where B z0 is the typical reconnected 
magnetic field. In summary, in the relativistic regime, the AR (Fig. [l2|c) becomes larger 
due to its relativistic outflow speed v out /c, and then reconnection generates more nonthermal 
particles (Fig. [l~2~|5). 



The direct acceleration around the AR is not likely to occur in the solar flare regions or 
the Earth's magnetotail, where the reconnection outflow is on the order of the Alfven speed 
of 10 6 (m/s) Cc~3x 10 8 (m/s) and the AR is restricted within the electron inertial length. 



4. Relativistic drift kink instability 
4.1. Simulation result 



Next, we introduce the other two-dimensional simulations in the y-z plane. The system 
size is now set to 1 (x) x 256 (y) x 512 (z). We set no artificial electric fields to excite 
instabilities, since instabilities spontaneously arise from thermal noise. We set no background 
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plasmas. In the first simulation, we set T = mc . This run is run D3 in Table p 



Characteristic snapshots of the simulation are presented in Figure 13 The color contours 



shows the plasma density, which is normalized by the peak density of the Harris current sheet 
p = 27^0. Figure 14 presents the energy spectra at the corresponding times. Figure 15 



presents the electric structures of selected stages. The left three panels of Figure [15] are color 
contours of the electric field E y . The sign of E y is positive in light gray (yellow) regions, 



while it is negative in dark gray (blue) regions. The right panels of Figure 15 show the 



electric field E z at the corresponding time. White lines are contours of B x , whose structure 
is similar to the density profile. 

After tens of the light transit time, we observe a linear growth of a kink-type structure 
around t/r c = 46.0. This is due to the relativistic drift kink instability (RDKI) (Zenitani & 



Hoshino| 2005a), a relativistic extension of the drift kink instability (DKI), which is a long- 
wavelength, current-driven instability in a thin current sheet. In the top panels in Figure 
15 one can see a typical structure of the polarization electric fields of E y and E z at this 
stage. The E y plot in the top left panel shows two types of color regions along the current 
sheet: yellow means E y > or in the right direction and blue means E y < or in the left 
direction. They are anti-symmetric with the neutral plane (z = 0). Since B x is positive 
on the upside of the current sheet and negative on the downside of the current sheet, the 
observed sign of E x B is consistent with the ^-displacement of the plasma bulk motion. The 
other components of electromagnetic fields, E x , B y , and B z , are negligible. After t/r c ~ 50, 
the system turns into its nonlinear stage. The typical signature is presented at t/r c = 64.0 
in Figures 13 14, and [15} In this stage, the current sheet is strongly modulated by the 



RDKI around —4 < z/X < 4. Because of the ^-displacement by the RDKI, we see the 
characteristic regions in a row along neutral plane (z = 0) in the middle left panel of Figure 
15 Originally, these regions were alternately arranged on the upper side and on the lower 



side of the neutral plane in the top left panel of Figure 15 Importantly, we observe a clear 



sign of the particle acceleration in high-energy tails in the energy spectra in Figure [14J We 
discuss the acceleration process later. 

After the nonlinear development, the folded current sheets start to collapse each other, 
and then the system turns into the "mixed" stage. Typical signatures are presented in 
snapshots at t/r c = 82.0 in Figures 13, 14 and 15 Particle acceleration seems to stop 



because the acceleration channel structure disappears. Because of the collisions of sheet 
fragments, magnetic energy is diffused into plasma heat, and then the total kinetic energy 
has increased by 170% at this time. Finally, the mixed current sheet slowly evolves into 
the broadened current sheet, which is 3 or 4 times thicker than the initial sheet width. 
The signature of this stage is presented in the snapshots at t/r c = 200.0 in Figures 13 and 



- 15 - 



14} The spectrum looks unchanged from the mixed stage; one can still see the remnant of 



the nonthermal tail. We still recognize both a remnant of the nonthermal tail and heated 
plasma component in the spectrum at t/r = 200.0 in Figure 14 Throughout Run D3, the 
total energy is conserved within an error of 0.3%. 



4.2. Linear analysis 

Throughout the simulation run, we have observed the perturbed magnetic field 5B X in 
the neutral plane (z = 0). We have examined the periodic perturbation modes and their time 



histories are presented in Figure 16 Each mode (2-5) corresponds to k y X = 0.49,0.74,0.98, 
and 1.22, respectively. The most dominant mode is mode 4 with k y X = 0.98. Its linear 
growth rate is u>i ~ O.llr^ 1 or oji ~ 0.035fi c , where fl c = oj c /^ = (eB) / (■jmc 2 ) is the typical 
gyro frequency. 



Pritchett et al. (1996) (hereafter P96) analyzed the linear stability of the DKI with 
fluid equations and Maxwell equations. We employed a similar method to our case with the 
relativistic fluid equation 

^{p s + e a ){— + v s • V)v s = -Vp s + ls q s n s [E + -j x Bj - -f ( 7s g s n s E • v s + -^), (18) 
the particle conservation law 

d 

^(j s n s )+V-(j s n s v s ) = 0, (19) 

and adiabatic gas condition with polytropic index T = 4/3 

Ps oc n r s . (20) 

The polytropic index T should go down to 5/3 in the nonrelativistic limit; however, setting 
T = 5/3 provides no significant change in our calculation. The detailed formulation is pre- 
sented in Appendix B. We have to be careful because the typical scale of particle motions 
are comparable to the plasma sheet thickness: the gyroradius c/(eB/^mc) ~ 0.3A or the 
meandering width [(c/n c )A]~ 1 / 2 ~ 0.6A. Therefore, fluid theory and the adiabatic gas con- 
dition are no longer valid in small scale structures. In addition, the theory is less trustable 
when we choose a larger f3 parameter, because the typical gyro scale becomes larger (eq. 
[To"]). The eigenprofiles for k y \ = 1.0, which is close to the most dominant mode (mode 
4; k y X = 0.98) in run D3, are presented in Figure 17 In the top panel in Figure 17, per- 
turbations of the magnetic field (5B X ) and the electric fields {SE y and SE Z ) are shown as 
functions of z. They are normalized by the maximum amplitude of 5B X . The anti-symmetric 
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structure of 5E y corresponds to the positive and negative E y regions along the neutral plane 
(z = 0). The other three components (5B y ,5B z ,5E x ) are zero in this case. The ^-profiles 
of velocity perturbations are presented in the middle panel in Figure 17 The parameters 
SV y ± = 5v yp ±5v ye and SV Z ± = Sv zp ±5v ze show the bulk or relative motion of the two species. 
The x components of the velocity perturbations {8v x s pfi \) are zero. The plasma bulk motion 
in the z-direction, 8V Z+ = Sv zp + Sv ze , is dominant around the neutral plane. The relative 
^-motion V z _ = 5v zp — 5v ze is also dominant around the neutral plane. It is mainly due 
to z-projection of the original counter-streaming current. A small-scale structure around 
\z/X\ ~ 2 — 3 may be numerical noise, which we could not completely clear out. The bottom 
panel in Figure 17 shows the density perturbation of 5D± = 5d p ± 8d e , where d = jn is the 
plasma density in the observation frame. The anti-symmetric structure of 5D + is consistent 
with the z-displacement of the central current sheet. Generally, their structures are in good 
agreement with the obtained simulation data, except for unphysical small structures (i.e. 
small waves in i5V y + around z/X = ±2.5). The perturbation profiles in the simulation at 



t/r c = 40.0 and the eigen functions are presented in Figure 18 The left panel compares 
mode 2 (kX = 0.49) in simulation and eigen mode for kX = 0.50. Both are in excellent 
agreement, although SE y in simulation looks slightly noisy because of its weak amplitude. 
The right panel of Figure 18 compares mode 4 (kX = 0.98) and eigen mode for kX = 1.00. 
They are in good agreement in SB X , while 5E y in simulation is approximately half of the 
eigen profile. The difference is also found in their the growth rates; the eigen growth rate at 
kX ~ 1.0 is uJi = 0.19T" 1 . We believe the fluid instability is dumped by the kinetic effect, or, 
in other words, fluid theory is losing its validity around the short wavelength of kX > 1. We 
discuss the validity of our theory later again. The charge separation (D-) is schematically 
described in Figure 19 Once the current sheet is modulated, the streaming e + fluid and the 
counter-streaming e~ fluid are slightly separated in ±x direction, due to the particle inertia 
or the drift motion between the centrifugal force and magnetic fields. Then the separated 
charges produce the electrostatic fields E y and E z . Since the E y x B x motion enhances the 
current sheet modulation in the z direction, the instability continues to grow. There is no 
influence of the Kelvin-Helmholtz instability in the system. Since we set empty background 
plasmas, there is no velocity shear between the current sheet and the background region. We 
also compare run D3 and the other run with 5% of background plasmas (run D3a in Table 
|3|), and we confirm that their growth rates are the same. 



4.3. Analytical theory 



Following P96, we construct a simple analytic solution of the instability in a long- 
wavelength limit and with ^ < 1. We drop the third term from the right-hand side of 
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equation (18), because they are smaller than other terms when linearized in a limit of /? <C 1. 
We get 



j'Vs/c 2 )('y s q s n s SE ■ v 5 
%q s n s 5E 



Q(n (v s /c 2 )d5p s /dt ^ < 



0((3 . 
c 



A) 



(21) 



where c s is the sound speed. In addition, using the density (n s , y s = d s ) in the observer's 
frame and rewriting the enthalpy term M = m r y s {l + [T/ (T — l)][T/(mc 2 )}), we obtain 



dt 



Md s (— + v s ■ V)v s = -X7 Ps + d s q s (e + — x B 

V c 



and the equation of continuity 



dds 
~dt 



+ V • (d a v) = 0. 



(22) 



(23) 



Using P96's method (described in the Appendix in P96) in equations (22) and (23) and 
the Maxwell equations, one can obtain the modified growth rate of the DKI, 



lRDKl/(k y Cp) » 1 -2^ 



(eB /Mc)X 



k y X = l- 2 7/3 (l + 



T 



c(3 



T — lmc 2 \(eBo/mc)X 



k y X, (24) 



where 'Jrdki is the growth rate of the RDKI. Since equation (24) gives larger results than 
the actual growth rates in nonrelativistic studies (Pritchett et al. 1996; Daughton 1999), 
we assume that its right-hand side of the equation gives an upper limit of the growth rate. 
Using equations (11) and (12), we obtain Jrdki as a function of k y X, 



/ Trie 2 r \ 
TcIrdki « [1 - lP\-^r + p— J£ 2 (M)] • 



(25) 



This equation has a maximum value 



Tel RDKI 



1 /mc 2 r 



-i f (4 7/3 /3)- 1 (T/mc 2 ) (T < mc 2 ) 
(167/3/3)- 1 (T > mc 2 ) 



(26) 



at 



kyX 



mc r 



r - 1 



-1 



(2 lf3 (3 2 )- 1 (T/mc 2 ) (T < mc 2 ) 
(8 7 ^ 2 )- 1 (T > mc 2 ) 



(27) 



Thus, as long as k y X < 1 is satisfied in equation (27), equation (26) gives a plausible upper 
limit of Jrdki- In addition, equations (p6j) and (27) can be rewritten as 



t, 



Hrdki < (kX y /2)p < 13. 



(2? 
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When k y X > 1 in equation (27), we use k v X = 1 in equation (25) to obtain a plausible 



extension of the upper limit, considering that the original equation (24) is no longer valid in 



k y X > 1. We have 



/mc 



TdRDKI = [1 - 1P{-JT + Y^l)^ ] ' 13 



[1 - -ip(3 2 (mc 2 /T)}(3 (T < mc 2 



A lp (3 2 )(3 



(T > mc 2 ) 



Note that both equations (28) and (29) satisfy the simple condition 

Tel RDKI < P- 



(29) 



(30) 



Figure 20 compares the linear growth rate of the instability in units of the light transit 
time of t c = X/c (left axis) or the typical gyro frequency Q c (right axis), as a function of the 
normalized wavenumber (k y X). The calculated eigen growth rates, the analytical estimate 
(eq. [25]) and the observed growth rates in the simulations are presented. The observed 



growth rates correspond to modes 2-5 in Figure 16 We discuss the RDSI line in §4.4. The 
three values are in good agreement when k y X < 0.7. In shorter wavelengths of k y X > 1, 
the simulation results do not agree with the theories because the instability is suppressed 



by kinetic effects. Daughton (1999) investigated the nonrelativistic DKI using Vlasov code, 



and his result also agreed with the fluid theory, when k y X < 0.7. 



4.4. Relativistic drift sausage instability 



The RDSI line in Figure 20 stands for the relativistic extension of the drift sausage 



instability, a cousin mode of the drift kink instability (Biichner & Kuska 1999 Yoon & Lui 



2001 |Silin &: Biichner 2003). Hereafter, we call it the relativistic drift sausage instability 
(RDSI). Its eigen profiles at kX = 0.25 are presented in Figure 21 The anti-symmetric profile 



of 5B X and the symmetric profile of 5D + are consistent with its sausage-type structure. The 
RDSI is usually faster than the RDKI in a very long wavelength limit of kX <C 1. On the 
other hand, the fastest RDKI is usually faster than the RDSI. In our runs, we could not 
observe any sign of the RDSI. This result is consistent with P96, who suggested that Landau 
damping stabilizes the (nonrelativistic) drift sausage instability. 



4.5. Particle acceleration 



In this section, we discuss the particle acceleration in the nonlinear stage, based on 
our previous works (Zenitani & Hoshino 2005a). To study the acceleration site, we select 
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high-energy positrons whose y-positions (y) and kinetic energies (e) at t/r c = 82.0 satisfy 
the conditions 6.4 < y/X < 12.8 and e/mc 2 > 20. Their spatial distributions and energy 
spectrum at three stages (light gray at t/r c = 46.0, gray at t/r c = 64.0, and dark gray 



at t/r c = 82.0) are presented in Figure 22 The selected positrons are accelerated in an 



'acceleration channel" (AC) around z ~ into the +x direction. Inside the AC, positive 



E y regions are in a row (e.g. middle left panel in Fig. 15), and then high-energy particles 
successively skip across the multiple positive E y regions along the AC, once their Larmor 
radii are comparable with or larger than the quarter wavelength of the RDKI. In energy 



spectra at t/r c = 46 and 64 in Figure 14, particle acceleration effect is only observed in the 
high-energy range of e > 10 — 12mc 2 , where the gyro radius (ri = -y(c/uj c ) > (1.5 — 2)A) 
exceeds the quarter wavelength n/(2k y ) ~ 1.6A. This threshold condition "fc/u c > n/(2k y ) 
yields 

e > 7reB \/(2k y \) = ttT/((3 k y X) ~ ttT//3. (31) 

So, the particle acceleration in the AC will smear out in the thick current sheet of f3 <C 1, 
because only a few particles satisfy this criterion. 

The folded structure evolves into the mixed stage when its ^-displacement (Az) is nearly 



the same as the half- wavelength of the RDKI (Az ~ n/ky ~ ttA). So, in Zenitani & Hoshino 



(2005a), we evaluated the timescale of the nonlinear stage (tjv) as 

t n ~ n\/v z ~ (tc/P)t c , (32) 

where we assume that the typical z- displacement speed v z ~ (0.2 — 0.3)c is approximated 
by the typical counter-streaming speed of [3c. After that, a lot of low-energy particles in the 
current sheet start to interact with the electric fields. The magnetic energy consumption 
is well explained by Joule heating between the average electric field E y and the zeroth 
order current J y ~ cB /(4:iiAz) in the broadened current sheet. We also approximated the 
dissipation timescale (td) considering the energy consumption rate: 

r D ~ {2AzBl/%u){2AzEyJ y )- 1 ~ (n/(3)r c , (33) 

where E y ~ (1/2)(v z /c)Bq and where the factor of (1/2) represents that the electric field 



is not uniform around the AC. From equations (11), (12), (|32[) , and (33), we estimate the 



maximum energy gain during the nonlinear stage and the mixed stage, 

Ae est ~ ecE y {r N + r D ) ~ eB^X ~ 2ttT//3. (34) 



Typical estimated values e est = e max0 + Ae est are presented in Table [3j 
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4.6. Case studies 

We have further carried out simulation runs with various temperatures, runs D1-D5 
(T/mc 2 = 16,4, 1, 1/4, and 1/16) in Table [3] We observe current sheet modulation and the 
electromagnetic perturbations by the RDKI in all five cases. Relativistic runs (T/mc 2 > 1) 
show similar results in many aspects of energy spectra, wavelength, and growth rates. The 
energy spectra in runs Dl and D2 have both the hot thermal component and the nonthermal 
component and they look "proportional" to those in run D3. The wavelength and the growth 
rates are nearly the same, when normalized by the sheet thickness or the light transit time. In 
less relativistic runs (T/mc 2 < 1), we observe longer wavelength modes (mode 2; k y X = 0.49) 
as predicted by the theory (eq. [27]). A high-energy tail becomes less apparent (Figure [23^) 
and the RDKI grows slower than the theory (eq. [26]). This is probably because the fluid 
theory loses validity in the nonrelativistic limit (eq. [l~6]). 



Figure [23] a shows enhancements of the particle kinetic energies or the released magnetic 
energy in runs D1-D5 with T/mc 2 = 16,4,1,1/4,1/16. Energies are normalized by the 
pressure energy in the original current sheet. The saturation level always seems to be 5 — 6, 
and the timescale for the energy release is t ~ 20r c . Nonrelativistic cases take more time due 



to their slower growth rates. Comparing with reconnection cases in Figure 12 s, the amount 



of released energy is comparable; however, due to its faster growth rates their timescales 



are shorter in the relativistic regime. Figure [23p presents the nonthermal ratio parameter of 
K n th/K. We could not resolve the nonthermal part with our method in run D5 (T/mc 2 = 
1/16) and in the late-thermalized stage of run D4 (T/mc 2 = 1/4). The parameter seems 
to become small as a smaller T is used and the particle acceleration is enhanced in the 



relativistic regime of T/mc 2 > 1. Comparing Figures 12b and 23!*, apparently reconnection 
produces more nonthermal energy than the RDKI. Although the width of the simulation box 
is not the same (25.6A[y] for the RDKI, 102.4A[x] for reconnection), the RDKI's nonthermal 
ratio will not change. If we look at the smaller region around the X-type region, we will 
obtain a higher nonthermal ratio. Remember the partial spectrum at t/r c = 100.0 in Figure[3] 
Apparently, the partial spectra is highly dominated by nonthermal energy and it is integrated 
over —12.8 < x/X < 12.8 — the same width as the RDKI simulation domains (25. 6A). So, 
as an origin of the nonthermal particles, reconnection is a more favorable candidate due to 
the acceleration processes around the X-type region. On the contrary, the RDKI quickly 
converts magnetic energy into plasma thermal energy or plasma heat. 

The maximum energies (e max ) obtained by the RDKI are presented in Table [3] They are 



compared with our estimate from equation (34), which seems to be a good approximation 



for these runs with (3 = 0.3. This relation will not be valid in the "thick" current sheet with 
small (3. In thicker current sheets, the Larmor radius becomes smaller compared with the 
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wavelength, so that particle acceleration may not take place. The maximum energies (e max ) 
obtained by reconnection are also presented in Table [3] As far as we have investigated, 
maximum energies by reconnection are far larger than those by the RDKI, and reconnection 
is undoubtedly favorable to accelerate high-energy particles. By the way, we do not know the 
upper limit energy of acceleration by the single reconnection. In Table [3j we roughly estimate 
the maximum energy by e est = eE y c(L x /2c) = eB L x /2 considering that the acceleration 
electric field E y is on the same order as B Q and considering that reconnection continues 
until two outflow jets toward the indirections meet each other at the periodic boundary at 
x = ±L x /2. 



5. Discussions 
5.1. Parameter dependence and comparison 



Figure 24 shows theoretical and observed growth rates of simulation runs R1-R5 and 
D1-D5 with different temperature parameters T/mc 2 = 1/16, 1/4, 1,4, 16. The growth rates 
(uji) are normalized by the light transit time r c . The linear growth rates of magnetic recon- 
nection are calculated from the reconnected field energy B 2 in the neutral plane. The 
normalized growth rates are around T c u>i ~ 0.03. Theoretical growth rates of reconnection 



are represented by the theory of the relativistic tearing mode (Zelenyi & Krasnoselskikh 



1979). Its growth rates 'Jrti in pair plasmas are: 



4= k\(l - k 2 \ 2 )pl 2 ^Thpmc 2 )- 1 ' 2 (T«mc 2 ) 
1 2 -f k\(l - k 2 \ 2 )pl 2 (T > mc 2 ) { ] 

In the relativistic case, we obtain the maximum growth rate by setting k\ = l/\/3 

TdBTI < 0.35/? 3 / 2 , (36) 

For (3 = 0.3, the maximum growth rate is t^rti ~ 0.055. This rate is 1.7 — 2.0 times faster 
than simulation results, but remember that reconnection is not identical to the tearing mode. 
The linear growth rates of RDKI/DKI are obtained from the perturbed magnetic field 5B X 
in the neutral plane, and the growth rates of the most dominant mode are selected. The 



growth rates of the RDKI/RDSI are estimated by equations (26) and (29). The growth rate 
of the DKI increases as T increases in the nonrelativistic regime of T/mc 2 1. On the 
contrary, the RDKI is insensitive to T. The theory (eq. [26]) clearly explains this signature. 
Relativistic pressure enhances effective inertia through the fluid enthalpy term, and then 
the enhanced inertia cancels further growth of the RDKI. So, the rate stays constant in the 
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relativistic limit. Since equation (27) slightly exceeds unity (k y X ~ 1.1 — 1.3) in relativistic 
cases of T > mc 2 , another criteria of equation (29) is presented as a second branch of the 
theoretical rate in Figure 24 Both criteria give Tc'Jrdki ~ 0.19, which limits the obtained 
value of T r uj; ~ 0.11. 



Figure 25 shows eigen and theoretical growth rates as functions of the current sheet 
thickness (3 = Xp/X (eq. |1~4~]). In this case, the plasma temperature is fixed to T = mc 2 . 
The upper limit of the RDKI (3 (eq. [30]), the upper limit of the relativistic tearing mode 
0.35/3 3 / 2 (eq. 
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eigen growth rates for the relativistic tearing mode, and eigen growth 
rates for the RDKI are presented. For eigen growth rates for the RDKI, we present two cases 
with kX = 0.5 and kX = 1.0 because we do not know what wavelength is the most dominant 
and because the DKI theories use the long- wavelength assumption (kX < 1). In the range 
of 0.03 < kX < 0.7, it seems that the eigen growth rates of the RTI and the RDKI/RDSI 
are well limited by their theoretical upper limit. It is still unclear how the instabilities grow 
in the thick current sheet and whether or not they are limited by the upper limits. This 
should be investigated by larger PIC simulations or relativistic two-fluid MHD simulations. 



Theoretically, equation (28) makes sense as an upper limit of the RDKI, because (3 also 
represents the current (eq. [12]) and because the RDKI is a current-driven instability. 



Figure 26 compares the theoretical growth rates as functions of two parameters: T/mc 



and (3. We use equations (26), (29), and (35). Generally, the RDKI/RDSI grows faster than 
reconnection. Most parts of the frontmost flat region (/3 1 and T ^> mc 2 ) rely on equation 
(28) or (29). It seems that the RDKI's rate goes down around (3 ~ 1, but our fluid theory is 



no longer valid here, because the current sheet is too thin (eq. [Mj ). Also, the fluid theory is 
not valid around the backward region of (3 > 0.1 and T/mc 2 <C 1, where the Larmor radius 
is comparable with or larger than the sheet thickness A (eq. [16]). Except for these invalid 
regions, we can conclude that the RDKI is more likely to occur in a relativistic hot condition 
of T/mc 2 > 1. 

Importantly, in the relativistic regime of T > mc 2 our simulations show no drastic 
change. The various aspects are similar: the evolution of the reconnection or the RDKI 
structure, time history of total energy, nonthermal ratio, relative height of the nonthermal 
slope, and the maximum energy of the nonthermal tail. These results are quite reasonable 
considering the following facts. The characteristic speed is around the light speed ~ c (Table 
[T| and the typical Larmor radius relative to the sheet width is nearly constant (eq. [l5]). It 
is straightforward that the typical gyro period (fi" 1 = 7/^ c ) is constant when normalized 
by the light transit time r c = A/c. So, we can conclude that the physical processes are 
"similar" in these relativistic CclSGS, clS long as we consider kinetics. Therefore, by employing 
the moderately-relativistic simulation (T « mc 2 ), we may figure out ultra-relativistic kinetic 
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evolution, which is computationally more expensive; one has to resolve the gyro motion of 
the lowest energy particle (uj^ 1 ) while the system evolution is governed by the typical gyro 
period (fi^ 1 = j/u c ). If we start to consider radiation effects, we need another approach 
because the synchrotron radiation highly depends on the ultra relativistic Lorentz factor: 7 
or T/mc 2 . 



5.2. Pulsar wind problem 



Let us discuss a possible application to the striped pulsar wind ( Lyubarsky fc Kirk||2001 



Kirk & Skjaeraasen 2003). Our comparison shows that the RDKI is the most likely process 

2\ 



in the current sheet, as long as plasmas are relativistic hot (T 



> 



mc 



If the current-aligned 



magnetic field (the so-called guide field, B y ) exists, secondary magnetic reconnection will 



dominate (Zenitani & Hoshino 2005b); however, magnetic fields are expected to be highly 



toroidal inside the striped pulsar wind due to the fast rotation of the neutron star, and 
therefore the guide field would be weak or unlikely to exist. Therefore, here we assume that 
the RDKI dominates in the striped current sheets. 



Kirk & Skjaeraasen (2003) (hereafter KS03) has examined the magnetic dissipation prob- 



lem by using the expanding current sheets model. They estimated the magnetic dissipation 
rate, assuming that striped current sheets are expanding with a "dissipation speed" of c/3 c , 
where (3 C is the dimensionless speed. They assumed magnetic reconnection as the most likely 
process, and they estimated the dissipation speed with the growth rate of the relativistic 
tearing mode ( Zelenyi fc Krasnoselskikh||1979 ) in their "tearing mode-limited" scenario. Im- 
portantly, one- dimensional simulation results are well described by an asymptotic analytic 
solution in all cases. From the flux-freezing equation (eq. [38] in KS03), the entropy equa- 
tion (eq. [39] in KS03) and the dissipation speed (/? c ), the asymptotic solution A oc (r/ri)~ q 
can be obtained, where A is a ratio of the current sheet thickness to the wavelength of the 
striped current sheets and q is the spectral index. Two of their scenarios are presented in 
Table [2j They assume that the current sheets are completely dissipated at the radial distance 
of T max jr^c from the central neutron star, where ric is the distance of the light cylinder. 
Two columns of r max /ric are presented; the first column is the analytical form, and the 
second column is the estimated value by using Crab parameters. The termination shock is 
at r/r LC = 2.0 x 10 9 , and so r max /r LC < 2.0 x 10 9 is a favorable result. The "fast" scenario is 
a physical upper limit; its dissipation speed is based on the relativistic sound speed of c/v3. 
In the tearing mode-limited scenario, they set c(3 c = Xjrti, where t c ^rti = P 3 ^ 2 (eq. [35] 
in KS03) is a simplified form of equation (J35j) . If we modify the dissipation speed with the 
maximum growth rate of t c Jrti ~ 0.35/5 3 / 2 (from eq. 36 ), the dissipation distance (r max ) 
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increases by 50%. The modified value is presented as the "tearing mode-limited (modified)" 
scenario in Table [2] 

We propose a similar dissipation scenario by the RDKI. We consider that the current 
sheet evolves to the broadened current sheet by the RDKI and then it continues to expand 
by the cascading RDKI processes. We call it the "drift kink mode-limited" scenario. We 
assume that the current sheet expands by h ~(3-4) times by single RDKI process: the 
sheet's typical width A becomes A' ~ hX, where the prime denotes the physical value in the 
broadened current sheet after the first instability is saturated. Consequently, the plasma 
density in the current sheet is n' ~ tlq/H. For simplicity, we ignore the particle acceleration 



effect, because nonthermal particles carry less than 10% of kinetic energies (in Fig. 23?)). 
Assuming the pressure balance condition 2n T = 2n' T' = Bq/8tt, the second-stage plasma 
temperature increases to X" ~ hT due to the nonlinear mixing of the RDKI. Then, the 
Debye length (X' D /X') ~ (A^/A) and the Larmor radius (r' L /X') ~ {tl/X) remain constant, 
relative to the sheet width. Regarding the the relativistic Harris parameters, /3' ~ (3 remains 



constant, while the temperature increases, T' ~ hT. According to Figure 26, the RDKI is 
likely to occur again in the broadened condition. So, it is quite plausible that the current 
sheet is dissipated by series of the RDKI processes. Topologically, the RDKI is favorable 
for the current sheet dissipation scenario, because it does evolve into the broadened current 
sheet structure, while it is not sure that reconnection evolves into the broaden current sheet. 



In thick current sheets of (3 <C 1, we evaluate the dissipation speed using equation (30); 



c(3 c ~ X^rdki = c(3. It is a good approximation of the eigen growth rates by a factor of 1-3 



(Fig. 25). Using this, following KS03's method, one can obtain the asymptotic index q = 2/5 
and relevant values. These values are presented in Table [2j If we use the typical parameters, 
the obtained dissipation distance is r max /r L = 3.1 x 10 9 . This value is 3 times better than 
the original tearing mode scenario and the drift-kink-mode scenario substantially improves 
the dissipation model in the striped pulsar wind. 

By the way, our DKI/RDKI theory is based on a long- wavelength assumption. In a 
nonrelativistic ion-electron plasma, it is reported that short-wavelength current sheet insta- 



bility grows faster than long- wavelength modes (Suzuki et al. 2002). Magnetic dissipation 
becomes even faster, if short-wavelength mode grows faster in a thick current sheet. In a 
thin current sheet, the nonthermal accelerations are more likely to transfer magnetic energy 
into the y-momentum ±p y of nonthermal plasmas. Then, since the ^-pressure relatively 



decreases, the current sheet may be thinner and the RDKI may grow faster (Fig. 26). In 
a thick current sheet of /? <C 1, we can probably ignore nonthermal acceleration because 
few particles satisfy the threshold condition (eq. (3l]). When neighboring current sheets 
come close to each other, they start to collide each other during the nonlinear stage of the 
RDKI. We compare several runs with different simulation box size L z (runs D3, D3b, and 



- 25 - 



D3c in Table [3J). In run D3c, the current sheets collide with each other across the periodic z- 
boundaries and more than 90% of the magnetic energy is dissipated into particle energy. On 
the other hand, the actual dissipation speed may be somewhat slower than the growth rates 
of the instabilities (both the RDKI and the tearing mode). So, an effective dissipation speed 
should be further investigated by a larger simulation, in which we can observe cascading 
evolution of the RDKI. If the current sheet structure survives at the termination shock, the 
magnetic fields may further be dissipated by the collision of the shock and current sheets. In 
the downstream of the shock, magnetic reconnection may be triggered by turbulent electric 
fields. 

In summary, we propose the RDKI as an alternative magnetic dissipation scenario in 
the striped current sheets, based on our theoretical estimates. At present we do not know 
whether or not the RDKI is the final key to explain the a problem; however, the drift kink 
mode scenario is obviously more favorable than the reconnection/tearing mode scenario. 



5.3. Miscellaneous discussions 



In ion-electron reconnection, various instabilities occur in different scales. For exam- 
ple, the LHDI occurs on an electron scale, the DKI on an ion scale, reconnection is on a 
macroscopic MHD scale, and then coupling across these scales leads a complexity of plasma 
phenomena. In this study, we investigate a very thin current sheet in e plasmas. Conse- 
quently, both the reconnection and the RDKI occur in the same scale as the current sheet. 
Both of them greatly violate the current sheet structure, and we discuss them as compet- 
ing processes. When we consider a "thick" current sheet, the RDKI or other cross-field 
instabilities may occur on a rather microscopic scale, and then they may interact with the 
macroscopic reconnection process. These situations should be further investigated. 

It has been believed that the reconnection's global energy dissipation is controlled by 



the Hall physics (Birn et al. 2001). In pair plasmas, the relevant positron's contribution is 



canceled out by electron's, and so the Hall effect does not exist. However, simulations of pair 



plasma reconnection ( 


Zenitani & Hoshino 


2001 


Jaroschek et al. 


2004 


2005b Bessho & Bhattacharjee||2005 


) reported fast reconnection rates 



>0), 



Ey/ '(vaB )] of 0.1 or above. This issue is equivalent to the question "What is the origin 
of the reconnection electric field?". In classical nonrelativistic studies, the electric field is 
described by 



E 



-v e x B 



me 
e 



+ v e ■ Vv f 



(37) 



n e e e \ dt 

where v e is the electron's velocity and U e is the electron's pressure tensor. The Hall effect 
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comes from the difference in the electron Lorentz term (v e x B) and its ion counterpart, and 
it smears out in the case of pair plasmas. Bessho & Bhattacharjee (2005) showed that the 
pressure-tensor term explains the electric field in the vicinity of the neutral region. Similar 



analysis can be applied to relativistic cases, using the relativistic pressure tensor (Write & 



Hadleyfl975 ). In addition, the electron's inertia may be of importance due to the increasing 



mass from the Lorentz effect. 



5.4. Summary and Conclusion 

Let us summarize and conclude this paper. We have carried out series of PIC simulations 
of relativistic current sheet problems in pair plasmas. We examine the relativistic magnetic 
reconnection and we report that particles can be accelerated near the magnetic pile-up 
region as well as the X-type region. We have also studied the current sheet instabilities, 
which exists in the perpendicular plane to reconnection. The parameter survey shows that 
properties of these two processes are similar in the relativistic regime of T > mc 2 , as long 
as we consider the kinetics. In addition, by comparing their growth rates, the RDKI is the 
more likely process in the relativistic regime of T > mc 2 . Therefore we propose that the 
magnetic dissipation in the striped pulsar wind should be re-considered by using the RDKI. 
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A. Nonthermal ratio parameter 

In this appendix, we describe our "nonthermal ratio" index, to discuss the nonthermal 
acceleration. By integrating over half of the simulation box, we obtain the plasma energy 
spectra F sim (e). The total particle number iV and the total kinetic energy K satisfy 

poo 

N = / deF sim (e) (Al) 



o 



oc 



K = / de F sim (e)(e-mc 2 ). (A2) 
J o 

Next, we consider an ideal thermal gas, which is described by the Jiitter-Synge distribution 
with unknown temperature T. Similarly, we consider the particle number N gas , the total 
kinetic energy K gas and its energy spectra F gas (e) for this ideal gas. The Jiitter-Synge gas 



has the relation (Synge 1957) 



K gas /(N gas mc 2 ) = 3(T/mc 2 ) + [K 1 (mc 2 /T) / K 2 (mc 2 /T) - 1], (A3) 

where Ki )2 (x) is the modified Bessel function of the second kind. Then, we assume that 
Ngas = N and K gas = K, so that F gas represents an "equivalent" thermal distribution. 



Since equation (A3) is a monotonically increasing function of (T/mc 2 ), we can find a unique 
solution of T = T'. Then, we calculate the nonthermal parameter by comparing the original 
spectrum F sim (e) and the equivalent thermal spectrum F gas (e) with T = T'. The relation 
between the two spectra is illustrated in Figure 27 the spectrum F gas is presented in light 



gray, and the spectrum F S i m is presented behind it in white and dark gray. The two spectra 
cross at several points, e.g. e = £i,£2, and we define the nonthermal energy K nt h as 



K nth = I [F stm (e) - F gas {e)](e - mc 2 ) de, (A4) 

£2 



where F sim exceeds F gas at e = e 2 in the highest energy tail. This "nonthermal" part is 
presented in dark gray in Figure [27} We use the ratio of the kinetic energy carried by this 
tail (Knth) to the total kinetic energy (K) as the nonthermal ratio parameter. 



B. Linear analysis of relativistic drift kink/tearing mode 

We start with the relativistic fluid equation of motion, 

iSU = (^+p)«v+^ (bi) 
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where T(fi U id) is the energy-momentum tensor for fluid, e is the energy, p is the pressure, u 
is the mean four-velocity of the fluid, and r/ a/3 is the metric tensor. In order to preserve the 
energy and the momentum in the system, 



d 



(t; 



■a/9 



a/3 



Q x p\ (fluid) 1 (EM) 







(B2) 



should be satisfied, where T^em) is the electromagnetic stress-energy tensor. Using the 
relation 

fl^#0 = - (B3) 



where F is the electro-magnetic tensor and j is the four current, equation |B2| can be modified 
to 







■a/3 



Q x f3 (fluid) 



(B4) 



The spatial three components of equation (B4) can be described (See Sakai & Kawata (1980)) 
as 

,.2 



C 5 



d f v \ v dv 

{p s + e s ){— +v s ■ V)v s = -Vp s + %q s n s (E+ -± x Bj - -|( 7s g s ri s E ■ v s + -^), (B5) 



where the subscript s denotes the species, p for positrons and e for electrons. We also use 
the particle conservation law, 



d_ 
dt 



(j s n s ) + V • (7 s n s v s ) = 0, 



and Maxwell equations, 



„ „ 4ttv-^ IdE 

V x B = — y ~i s q s n s \ s H — 



19B 



(B6) 

(B7) 
(B8) 



VXE: 

We assume the plasma is an adiabatic gas, 

p s oc n T s . (B9) 

where T is the polytropic index of T = 4/3. Then, we replace some variables in equations 
p5b, pi, dB7b and p8b by using 



e s = n s m c + 



1 



Ps 



r- 1 

Is = [1 - (V S /C) 2 ]- 1/2 



(BIO) 

(Bll) 
(B12) 
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where d s stands for the number density in the observed frame. 



Finally, equations (B5), (B6), (B7), (B8) and (B9) are linearized for the the Harris 



field model, assuming that perturbations are given by 5f oc 5 f (z) exp(ik x x + ik y y + uit). 
The operators are transformed as d/dx — > ik x and d/dy — > ik y , and we normalize the 
first-order perturbations as B = B (F(z)x + b y y) + B b, E = B e, v s = ±c/3 + c/3v s , 
d s = d s0 f(z) + d s0 d s , and p s = p s0 f(z) + p s0 p s , where F(z) = tanh(>/A) and f(z) = 
cosh~ 2 (2;/A). The parameter b y = B y /B , the amplitude of the guide field, is set to 0. For 
numerical simplicity, some physical properties are calculated in the form of summaries or 
differences between the positron values and electron values, D± = d p ± d e , P± = p p ±p e , and 
V± = v p ± v e . We have 



Of 



dt bx 
di by 

ik) 



d_„ 
dt Cx 

a* 
di 6y 



d_ 
dt 



d_ 
Of 







'iD. 



d_ 

dt 

d_ 
dt 
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dV 

d_ 

dt 



V,- 



, . . d A 
-ck y [ie z ) -\- c e y 



d A 
dz 



ck x c y ck y & x 



d 



-A 1 f(z)V x - + ck y (ib z ) - Ct^&j, 

-i4ip+ + /(*) V) + - ckxWz 

-Aif(z)(iV x -) + ck y b x - ck x b y 



(B13) 
(B14) 
(B15) 
(B16) 
(B17) 
(B18) 



r lk„iiD- ) < lf( : ) | k x (iV x+ ) + k y (iV y+ ) + ^V z+ 



r ik„D + < }[( : j ( k x V x _ + kyVy- - ^(iV z . 



-c(3f(z)V z+ 

8k y b + + c 

-cf3f(z)(iV z „) 
d > 3 2 d ~ 



c/%\4_ + k x -f-P + + 2A 3 (tb z ) - A 3 b y (iV z . 



7(* 



(B19) 

(B20) 
(B21) 

(B22) 

(B23) 
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c dt 
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f{z)dz f(z) 
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(B25) 



(B26) 



(B27) 



(B28) 



where f'(z) = {d/dz)f(z) = -(2/A) tanh(z/A) cosh- 2 (z/A), A 1 = 4necf3d /B , A 2 = 
( Po /c/3)[(l}/c 2 )(n mc 2 + [r/(r - l)]^)]" 1 , and A 3 = (ed B o/ 'c)[( 7 |/ 'c 2 )(n mc 2 + [T/(T - 
l)]po)] -1 - Note that equations (B21), (B22), (B25) and (B26) still contain time derivative 
terms. We reorganize the six equations ( gig ), flB20j ), (|B2lJ, p22j ), ( [B25] , and fB26j ) to 
obtain (d/dt)P + ,(d/dt)(iP_), (d/dt)(iV v+ ), and (d/dt)V v -. 




1-IM 3 /? 3 7 |/c V -A 3 L3/[cf(z)\ 
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(B29) 
(B30) 



Finally, after replacing <9/<9t — > u>i, all the above relations can be solved as an eigen value 
problem of the following matrix. Solving this problem, we obtain the growth rate of the 
instability uji and the z-structure of the perturbations for given k x and k y . The matrix is 
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(B31) 



We set 241 grids in the z-direction and 20 grids per A; therefore, we solve — 6A < z < 
6 A. The size of the matrix is 3856 2 . The conductive- wall conditions are assumed at the 
boundaries. For selected cases, we use 361 grids in z, 5776 2 grids per a matrix and we 
checked the effect of the grid resolution and the boundary positions. In the case of (3 = 0.3, 
the grid number shows an error less than 1%. In the case of (3 < 0.1, we have to choose a 
high resolution configuration with 5776 2 grids. For numerical stability, we sometimes add a 
small amplitude of the (nonrelativistic) viscosity term to the fluid equation. The equivalent 
Reynolds number (R e ) is set to R e = p \cP/[i = 10 6 — 10 7 . The eigen modes for the tearing 
mode or the RDKI/RDSI mode are usually obtained without this viscosity term. We need a 
small amplitude of the viscosity term to solve short-wavelength RDKI modes (k x — 0, k y > 1) 
and the oblique modes (k x ^ 0, k y ^ 0). 



Fig. 1. — Snapshots of relativistic reconnection at typical stages: t/r c = 60.0,80.0, and 
100.0. In the left panels, the solid lines represent magnetic field lines and the color contour 
shows the density of plasmas normalized by the original plasma sheet density p = I^^tiq. 
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Fig. 2.— Late-time snapshot at t/r c = 300.0 of run R3. Main (-12.8 < z/X < 12.8) and 
sub (12.8 < z/X < 38.4) simulation boxes are presented. The white rectangle shows the 
region of interest in Fig. [TJ 
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Fig. 3. — Particle energy spectra in run R3 are presented in log-linear and log-log format. 
The partial energy spectra at t/r c = 100.0 is also presented as t/r c = 100.0'. 
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Fig. 4.— Structure of the electric field in run R3. The left half displays the Lorentz 
invariant E 2 — B 2 , normalized by B 2 . The right half displays the reconnection electric field 
E y , normalized by B . 
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Fig. 5.— Spatial distributions of accelerated particles whose energy exceeds 50mc 2 at 
t/r c = 100.0. White lines show magnetic field lines. The red marks in the bottom panel 
show the spatial distribution of highly accelerated particles, whose energy exceeds 90mc 2 . 
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Fig. 6. — Particle trajectory (solid line) of a highly accelerated positron, its xz- and xy- 
projections (dashed lines) are presented. Labels in the xy plane show the relevant times. 
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Fig. 7. — Time history of the particle's energy, momentum and the electromagnetic field at 
the positron's position in the simulation frame. 
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Fig. 8. — x-t diagram of the positron trajectory and the reconnection electric field E y along 
the neutral sheet (z — 0). 
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Fig. 9. — Field structure in the neutral sheet at t/r c = 100.0. The normalized electromag- 
netic fields B z and E y and the ratio of |^|/|i? 2 | are presented. 




Fig. 10. — Snapshots of runs for T/mc 2 = 16 (top) and T/mc 2 = 1/16 (bottom) 



-44- 




Fig. 11- 



Particle energy spectra in run Rl in log- log format. 
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Run 


T/mc 2 


Vin_max / C 


Vout_max/ C 


Va/c 


yinflo W/c 


Rl 


16 


0.6 


0.97 


0.57 


0.99 


R2 


4 


0.4 


0.9 


0.57 


0.99 


R3 


1 


0.39 


0.88 


0.53 


0.98 


R4 


1/4 


0.29 


0.63 


0.45 


0.94 


R5 


1/16 


0.19 


0.53 


0.30 


0.84 



Table 1: Maximum outflow velocity compared with two types of the Alfven velocity: (1) 
the typical Alfven velocity and (2) the inflow Alfven velocity. 
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Fig. 12. — (a) Enhancements of the particle kinetic energies in the reconnection runs with 
T/mc 2 = 16, 4, 1, 1/4, 1/16. Energies and time are normalized by the pressure in the original 
current sheet and the light transit time, respectively, (b) Nonthermal ratio parameters for 
reconnection. (c) Size of the acceleration region. 




Fig. 13.— Snapshots of run D3 in the yz plane at various stages (t/r c = 
0.0, 46.0, 64.0, 82.0, 110.0, and 200.0). 




Fig. 14. — Energy spectra in run D3 at various stages. The particle count number, integrated 
over the main simulation box, is presented as a function of the particle energy. 




Fig. 15. — Electric fields E y and E z at t/r c = 
of B x . 



46.0, 64.0, and 82.0. White lines show contours 
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Fig. 16. — Time histories of the amplitude of perturbed magnetic field 5B X in the neutral 
plane [z = 0) for drift kink modes. Modes 2-5 correspond to k y \ = 0.49, 0.74, 0.98, and 1.22. 
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Fig. 17. — The ^-profiles of the eigen functions of the relativistic drift kink instability (RDKI) 



for the parameter k y \ = 1.0. Top: 



Perturbations of 5B X , 5E y and SE Z 



Middle: Perturbed 



plasma velocities. The bulk velocities of plasmas are 5V{ V:Z y + = 5v{ y:Z y p + 5v{ yjZ y e and the 
relative velocities of the two species are 5V{ ytZ y- = 5v{ y ^y p — 8v{ y , z } e . Bottom: The density 
perturbations are presented: 5D + = 5d p + 5d e and 5D_ = 5d p — Sd e . 
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Fig. 18. — Perturbation profiles of the simulation data and relevant eigen modes. Left: 
Mode 2 (kX = 0.49) and the eigenmode for kX = 0.50. Right: Mode 4 (kX = 0.98) and the 
eigenmode for kX = 1.00. 




Fig. 19- 



Schematic explanation of charge separation in the RDKI 
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Fig. 20. — Growth rates of the RDKI as a function of the normalized wavenumber 
{kyX). Eigen growth rates {solid line with plus symbols), another branch of eigen growth 
rates {dashed line with crosses), eq. ( [25] ) {dashed line), and simulation results {asterisks). 
The timescale is presented in units of r c {left axis) and fl c {right axis). 
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Fig. 21. — Same as Fig. 17 but for the relativistic drift sausage instability (RDSI) with the 
parameter k y \ = 0.25. 



Fig. 22. — Spatial distributions and energy spectra of selected positrons, which satisfies 
the following two criteria: (1) the kinetic energy satisfies e/mc 2 > 20 at t/r c =82.0 and (2) 
the ^/-position satisfies 6.4 < y/X < 12.8 at tjr c =82.0. Three stages of the distribution are 
presented, t/r c = 46.0 (light gray or yellow), 64.0 (gray or light orange) and 82.0 (dark gray 
or deep orange). 
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Fig. 23. — (a) Enhancements of the particle kinetic energies in the RDKI runs. Energies 
and time are normalized by the pressure in the original current sheet and the light transit 
time, respectively, (b) Ratio of nonthermal energy in the case of the RDKI runs. The ratio 
for T/mc 2 = 1/16 cannot be calculated by our method. 
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Fig. 24. — Temperature dependency of the growth rates (r c uji). Theoretical growth rates 
of RTI and RDKI, and observed growth rates of magnetic reconnection and the RDKI are 
presented. 
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Fig. 25. — Growth rates of the instabilities as a function of the current sheet thickness 
P = Ad/A. The growth rate u>i is presented in a timescale of r c . Upper limit of the RDKI j3 
(from eq. 30 ), upper limit of the relativistic tearing mode 0.35/? 3//2 (from eq. 35 ), obtained 
eigen growth rates for the tearing mode, obtained eigen growth rates for the RDKI and eigen 
growth rates for the RDKI for short wavelength are shown. 
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Fig. 26. — Theoretical estimate of the growth rates Ui, as a function of the plasma tempera- 
ture T/mc 2 and the drift speed (3. The dashed surface shows linear growth rates of magnetic 
reconnection in electron-positron plasmas, based on Zelenyi & Krasnoselskikh (1979). The 
solid surface shows the upper limit growth rates of the drift kink/sausage instability. 
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Table 2. Simple estimates of KS03 scenarios 



Dissipation 


q 


T*max J^LC 


fmax/fLC 


A(r/r LC )-i 


Tearing-mode (modified) 


5/12 


0.75 • M 4/5 i 3/10 


9.3 x 10 9 


[5 s 0.35 2 /(19 7 12tt 2 )] 1 /i2^-i/3| / -i/8 


Tearing-mode (KS03) 


5/12 


0.5 ■ M 4 / 5 L 3 / 10 


6.1 x 10 9 


[5 8 /(19 7 127r 2 )] 1 /iV- 1 /3L-i/8 


Drift-kink-mode 


2/5 


0.25 • v^L 1 / 4 


3.1 x 10 9 


4-3/8^.1/5^-2/5^-1/10 


Fast (KS03) 


1/3 


0.1/. 2 (l-/? 2 )//3 c 


4.6 x 10 7 


{6^/(25^(1 - pm i/3 » 2/3 



Note. — The dissipation index q, the dissipation distance r max , its estimated value, and the factor for A 
are shown. The Crab parameters fi = 2 X 10 4 and L = 1.5 X 10 22 are the total energy carried by the wind per 
unit rest mass (eq. [23] in KS03) and the dimensionless flow luminosity (eq. [29] in KS03), respectively. 
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Fig. 27. — Relation between the observed distribution F sim (white and dark gray) and the 
equivalent thermal distribution F gas [light gray). Two parameters E\ and 62 represent crossing 
points of the two spectra. The high-energy tail region in dark gray carries the "nonthermal" 
energy. 
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